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The particle-hole Density Matrix Renormalization Group (p-h DMRG) method is discussed as a 
possible new approach to large-scale nuclear shell-model calculations. Following a general descrip- 
tion of the method, we apply it to a class of problems involving many identical nucleons constrained 
to move in a single large j-shell and to interact via a pairing plus quadrupole interaction. A single- 
particle term that splits the shell into degenerate doublets is included so as to accommodate the 
physics of a Fermi surface in the problem. We apply the p-h DMRG method to this test problem 
for two j values, one for which the shell model can be solved exactly and one for which the size 
of the hamiltonian is much too large for exact treatment. In the former case, the method is able 
to reproduce the exact results for the ground state energy, the energies of low-lying excited states, 
and other observables with extreme precision. In the latter case, the results exhibit rapid expo- 
nential convergence, suggesting the great promise of this new methodology even for more realistic 
nuclear systems. We also compare the results of the test calculation with those from Hartree-Fock- 
Bogolyubov approximation and address several other questions about the p-h DMRG method of 
relevance to its usefulness when treating more realistic nuclear systems. 

PACS numbers: 21.60.Cs, 05.10.Cc 



I. INTRODUCTION 

The nuclear shell model Q| is one of the most powerful 
approaches for a microscopic description of nuclear prop- 
erties. The low-energy structure of a given nucleus is 
described in this approach by assuming an inert doubly- 
magic core and then diagonalizing the effective hamilto- 
nian within an active valence space consisting of at most 
a few major shells. Despite the enormous truncation in- 
herent in this approach, the shell-model method as just 
described can still only be applied in very limited nu- 
clear regimes, namely for those nuclei with a sufficiently 
small number of active nucleons or a relatively low de- 
generacy of the valence shells that are retained. The 
largest calculations that have been reported to date are 
for the binding energies of nuclei in the /p-shell through 
64 Zn §. 

For heavier nuclei or nuclei farther from closed shells, 
one is forced to make further truncations in order to re- 
duce the number of shell-model configurations to a man- 
ageable size. The most promising approach now in use 
is to truncate on the basis of Monte Carlo sampling || . 
In this way, it has recently proven possible to extend 
the shell model beyond the /p-shell to describe the tran- 
sition from spherical to deformed nuclei in the Barium 
isotopes j|. 

Another attractive possibility is provided by the Den- 
sity Matrix Renormalization Group (DMRG), a method 
that was initially developed and applied in the framework 
of low-dimensional quantum lattice systems H. For a 



recent review, see ref. Q. A simplified version of this 
method was applied M to the two-level pairing model, 
showing its convergence properties. Subsequently, the 
method was extended to finite Fermi systems and ap- 
plied to a pairing problem of relevance to ultrasmall su- 
perconducting grains H. This new methodology, which 
is referred to as the particle-hole (or p-h) DMRG, was 
recently applied to a first test problem of relevance to 
nuclear structure ||. The application involved identi- 
cal nucleons moving in a single large-j shell under the 
influence of a pairing plus quadrupole interaction with 
an additional single-particle energy term that split the 
shell into doubly-degenerate levels. Comparing with the 
results of exact diagonalization, it was shown that the 
method leads to extremely accurate results for the ground 
state and for low-lying excited states without ever requir- 
ing the diagonalization of very large matrices. Further- 
more, even when the problem was not amenable to exact 
solution, the method was seen to exhibit rapid exponen- 
tial convergence. All of this suggests that the DMRG 
method may indeed be a practical means of carrying out 
reliable large-scale shell-model calculations. 

Before proceeding to realistic applications of the 
method, it is important to clarify its various ingredi- 
ents and assumptions. Furthermore, it is critical to opti- 
mize how we apply the methodology for the subsequent 
more realistic and complex applications which will follow. 
With that in mind, we report here the results of a more 
thorough investigation of the recent application of the 
DMRG method to the one-orbit pairing plus quadrupole 
model. Several of the questions that were hinted at in 
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ref. H are now addressed, including the feasibility of ap- 
plying the method to very much larger model spaces. 

The outline of the paper is as follows. In Section II, we 
review the basic features of the DMRG method and pro- 
vide a fairly comprehensive description of the theoretical 
ingredients required for its implementation. In Section 
III, we describe the one-shell model that we will be using 
for the first tests of the DMRG method in nuclear physics. 
In Section IV we present the results of these test calcu- 
lations for a problem of 10 identical nucleon restricted to 
move in a single j = 25/2 shell. This problem can be 
solved exactly using the Lanczos algorithm, providing us 
with a good testing ground of the methodology including 
several variants not discussed in ref. In Section V, 
we discuss the application of the method to a problem 
involving 40 nucleons in a j = 99/2 orbit, a problem well 
beyond the limits of exact diagonalization. In Section 
VI, we summarize the principal conclusions of this work 
and outline some future directions for investigation. 



II. OVERVIEW OF THE P-H DMRG 
PROCEDURE 



A. Qualitative description 

The basic idea of the p-h DMRG method is to system- 
atically take into account the physics of all single-particle 
levels. The p-h nature of the method enters through 
a separation of the active single-particle levels into two 
sets that are separated by the Fermi surface; those above 
the Fermi surface are called the particle levels and those 
below are called the hole levels. The procedure begins 
by first taking into account the most important levels, 
namely those that are nearest to the Fermi surface, and 
then gradually including the others in subsequent iter- 
ations. At each step of the procedure, a truncation is 
implemented both in the space of particle states and in 
the space of hole states, so as to optimally take into ac- 
count the effect of the most important states for each of 
these two subspaces of the problem. The calculation is 
carried out as a function of the number of particle and 
hole states that are maintained after each iteration, with 
the assumption that these numbers are the same. This 
parameter, which we call p, is gradually increased and 
the results are plotted against it. Prior experience from 
other applications of the methodology suggests that the 
results converge exponentially with p. Thus, once we 
achieve changes with increasing p that are acceptably 
small we terminate the calculation. Assuming that the 
convergence is rapid, the method permits us to achieve 
an accurate description of the low-lying states of the sys- 
tem without ever having to diagonalize the large matrices 
that would arise in the absence of the DMRG truncation 
strategy. 



B. A bit more detail 



It is useful to now put the above qualitative remarks 
about the DMRG method into a more mathematical con- 
text. 

The general problem that we wish to solve is a stan- 
dard shell model problem of N v neutrons and N n protons 
(total N = N v + Nk) interacting via a one- plus two- 
body hamiltonian, with the neutrons and protons each 
restricted to an active set of single-particle orbitals. For 
simplicity of presentation, we will focus in this subsec- 
tion on just one type of particle, as is the case for the 
test applications we will report. In a later subsection, we 
will note those modifications required for neutron-proton 
systems. 

We begin by splitting the set of multiply-degenerate 
single-particle levels into an ordered set of doubly- 
degenerate levels. This is most naturally done by car- 
rying out an axially-symmetric deformed Hartree Fock 
(HF) calculation for the system. Each of the result- 
ing HF levels is then doubly-degenerate. Going from 
the set of spherical single-particle levels to the axially- 
symmetric set of doubly-degenerate levels is schemati- 
cally illustrated in figure 1 for a problem involving 8 
identical nucleons in the / — p shell. 



1/2 
5/2 
3/2 
3/2 

1/2 
7/2 



Pi 
f 



12 



5/2 
P 



3/2 



7/2 



1/2 
5/2 

3/2 
1/2 



FIG. 1. Schematic illustration of the splitting of the 
model-space single-particle levels into a set of dou- 
bly-degenerate levels by an axially- deformed Hartree Fock cal- 
culation. The dashed line represents the Fermi energy (Ef ), 
which separates the particle levels from the hole levels. Each 
doubly- degenerate level is labelled by its angular momentum 
projection on the intrinsic z-axis. 
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Each of these doubly-degenerate single-particle levels 
admits precisely four states. One is the state with no par- 
ticles, one is the state with one particle having positive m 
projection, one is the state with one particle having neg- 
ative m projection, and the last is the state containing 
two particles, with both positive and negative to. 

This set of doubly-degenerate single-particle levels nat- 
urally splits into two smaller sets, as likewise illustrated 
in figure 1. Those levels above the (non-interacting) 
Fermi surface are the particle levels and those below the 
Fermi surface are the hole levels. Reiterating what was 
stated earlier, we will systematically take into account all 
of the levels of the problem, by gradually moving away 
from the Fermi surface. 

In this basis, the shell model hamiltonian takes the 
form 

H ^ ^ e arn Ci^ arn CL arn 
am 

+ i (aiTOi, a 2 TO2|V r |a 3 TO3, a4TO 4 ) 

aim.ia 2 m 2 a 3 m3a^m^ 

X a a 1 m 1 a a 2 m 2 aa i m 4. a aarn 3 7 (1) 

with both a single-particle energy term and a two-body 
interaction. 

At this point, let us assume that we have already 
treated some of the levels of the problem, perhaps the 
first n particle levels above the Fermi energy (forming a 
particle block) and the first n hole levels below the Fermi 
energy (forming a hole block). Let us also assume that 
we have precisely p particle states in the particle block 
and p hole states in the hole block, i.e. the same number 
for each. 

Furthermore, we will assume that we have calculated 
and stored the matrix elements of all possible sub- 
operators of the hamiltonian within the p-dimensional 
particle and hole blocks, respectively. For a one- plus 
two-body hamiltonian, the operators of relevance are: , 
a^a, ala 1 , aWa and aJoJaa, as all others can be directly 
calculated from these using appropriate hermitian ad- 
joint relations. 

We now enlarge both the particle and hole blocks, 
by adding one more (doubly-degenerate) level to each. 
Clearly the number of states in the enlarged particle and 
hole blocks will be Ap. These states can be written as 
= where i = l,...,p and j = 1,...,4. 

What we would like to do is to truncate both the par- 
ticle and hole blocks from these 4p states to the optimum 
p states, precisely the same number that we had before 
adding the new level. 

How do we accomplish this? 

The first step is to recalculate the matrix elements 
of all of the hamiltonian sub-operators in the enlarged 
particle and hole blocks. To show how this is done, we 
will focus on one specific set of matrix elements, those 



of the operator a' a ro a a3mr Its matrix elements in the 
enlarged block can be readily expressed as 

(*> j\ a a 1 m 1 a a 2 m 2 |&j = (* l a iLimi a a 2 m2 I ^) Sj,l 

+ (jW aimi a a2 m 2 \l} 5i,k 

+ (-) n *<i|4 1 mJfc> UK 2 m 2 \l) 
-(-) nk {*K 2 m 2 \k) <j|< mi |0 , 

(2) 

where rn is the number of particles in state \i). Depend- 
ing on which of the two subspaces the operator indices 
refer (either the levels 1 — > n or the new level n + 1), a 
different term in the sum applies. For example, if both 
ati, tyi\ and 012, to 2 refer to single-particle states from the 
first n levels, it is the first term that applies. 

The key point is that all terms in (^|) involve a product 
of a matrix element from the previous p x p space (levels 
1 — > n) and a matrix element for the level that is be- 
ing added. The matrix element from the previous space 
was already calculated and stored, whereas the one from 
the new level is very simple and can be written down 
analytically. 

At this point we have calculated and stored all of 
the particle and hole matrix elements in their enlarged 
blocks. The next step is to construct the space of particle- 
hole states, by forming products of states in the particle 
and hole blocks, viz: \ip,ju) — \ip)\jh)- This is usually 
referred to as the superblock. 

We only build states in the superblock with the cor- 
rect number of total particles N for the problem of in- 
terest. That means that we only consider states with the 
same number of particles and holes. Furthermore, as- 
suming that the hamiltonian is rotationally invariant, it 
suffices to consider only those states for which the total 
z-projection of the angular momentum is zero. 

Within this product space, we then build the hamil- 
tonian matrix. Using the simplified notation 7 = {am}, 
this can be expressed as 

{ip,3H\H\k P ,l H ) =J^e 7 [{ip^a^kp) 8 jH j H 

7 

+ (iff|<4 a 7l^) ^P,k P ] 

+ i (liilj\ v \lk,li) 

7i,7j,7fc,7l 
Sip, k P 0ff| a 7, a 7 3 a 7i a 7 fc l^) 

+ 2{-) n "p {ip\al n \k P ) (jj?|4.a 7i a 7 J/ H ) 
+ 2(-)" fc p (ip\a 7l \k P ) (j H \a\.a\ jaih \l H ) 

+ (*p| a 7i a 7j l fc p) O'ff I 

~ 4(i P |a T4 a 71 |fcp) {j H \a\ o a lk \l H ) 

+ (ip|a 7i a 7 Jfc P ) (jffla^atjJjj) 

+ 2(-)" fc p {ip\a\.a 7k a 7l \k P ) (j H \a\.\l H ) 
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+ 2{-) n <*r(i P \a\ a* a 7! |fc P ) 0'*KM 



li 7j 7z 

(ip|a 7i a 7j a 7l a 7fe 



(3) 



Note that every term in eq. (|^) involves a product of a 
matrix element in the particle space and a matrix ele- 
ment of the conjugate operator in the hole space. All of 
these matrix elements, for both particles and holes, have 
already been calculated and stored in their respective en- 
larged blocks. 

Next we diagonalize the superblock hamiltonian, viz:, 



with 



l**> 



E *W\ip)\3B 



(4) 



(5) 



ip,jp = lAp 



Note that both sums go over the 4p states in the respec- 
tive enlarged particle and hole blocks. 

We now wish to construct the optimal approximation 
to the ground state of the system (or to some low-lying 
set of states) that is achieved when we only retain p par- 
ticle states and p hole states. By optimal, we will mean 
the approximation that maximizes the overlap between 
the truncated state and the exact ground state of the 
superblock in eq. (||)). For now, we will ask for 

the optimal description of the ground state; later we will 
discuss how to generalize this to several states. 

We will perform the optimized truncation in two 
stages, first asking what is the optimal approximation 
when we truncate the particle states and then asking 
what is the optimal approximation when we truncate 
the hole states. For simplicity of notation, we will let 
|*) = |*!) and *y = *g>. 

Let us focus first on the particle block. We wish to find 
a subset of particle states (p in number), such that the 
truncated ground state built up from these states has the 
largest possible overlap with the exact superblock ground 
state. To accomplish this, we first define the reduced 
density matrix for particles, 



p 

Pi pi' 



E 



(6) 



The reduced density matrix separates into blocks accord- 
ing to the number of particles and the m value. 
We then diagonalize this 4p x Ap matrix, 



> 



(7) 



A given eigenvalue co p represents the probability of find- 
ing the particle state \u P > in the full ground state wave 
function 1^ >. Thus, the optimal truncation in the sense 
described above corresponds to retaining the p eigenvec- 
tors | Up > that have the largest eigenvalues. 



We then do exactly the same for the hole block. 
Namely, we construct the reduced density matrix for 
holes, 



E 

ip = l,4p 



iff Vtl* , 

^IPJH ^ipj' 



(8) 



diagonalize it, and then keep only the p hole states with 
the largest eigenvalues. 

At this point, we have identified the p particle states 
and p hole states that best approximate the full ground 
state of the coupled superblock. The final step before 
proceeding to the next level is to transform all of the ma- 
trix elements that we calculated in the 4p-dimcnsional 
particle and hole spaces to the optimal p-dimensional 
truncated spaces. 

As noted earlier, we need not target our optimization 
to the ground state only. We can target it to any set of 
states we wish. If, for example, we wanted to optimally 
describe the lowest L eigenstates of the system, we would 
build mixed density matrices 

1 

I 



p 

Pi P i' r 



/ / L ' 1 p3h r p Jn 
k=l,L iff 



PjHj' f 



k=l 1 L ip 



(*0* 

ipjH ipj'. 



(9) 



and use them to choose the p most important particle 
and hole basis states to retain. 

Once this series of steps has been implemented, we 
simply return to the point where the discussion began 
and continue iteratively. Namely we add the next levels, 
precisely as we did above, first calculating and storing all 
matrix elements in the enlarged particle and hole blocks, 
then constructing and diagonalizing the hamiltonian ma- 
trix in the superblock, then constructing and diagonaliz- 
ing the reduced density matrices for particles and holes 
(depending on the states we wish to target) and then 
truncating to the most important p states for particles 
and holes, respectively, based on the eigenvalues of the 
associated density matrices. 

Were the number of particle and hole levels the same, 
we would simply proceed as above to treat all levels. We 
would initialize the iteration process by treating the two 
levels nearest to the Fermi surface, whose matrix ele- 
ments can be calculated very simply. We would then con- 
tinue to add levels as defined above. Truncation would 
not be necessary until the number of particlc/holc states 
exceeds p. Prior to that point, it is not necessary to con- 
struct the hamiltonian matrix in the superblock, since 
it is only needed for the purpose of truncation. Note 
that in the absence of truncation, the number of particle 
and hole states after treating n levels is 4". Thus, until 
4™ > p, truncation is not required. 

When there are different numbers of particle levels and 
hole levels, some modification to the above algorithm is 
needed. Assume for example that there are more particle 
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levels than hole levels, so that the system is less than half 
full. We can continue to add levels two at a time, one for 
particles and one for holes, until we have exhausted all of 
the hole levels. From that point on, however, there are 
no hole levels to add. Thus, in subsequent iterations, wc 
only add particle levels and only carry out the optimized 
truncation for those states. 

At each iteration, there are in fact 4p states for parti- 
cles and holes in the enlarged blocks required to build the 
superblock hamiltonian. It is only after truncation that 
it is reduced to p. Assuming that there are more particle 
levels than hole levels, we have a choice of how to proceed 
after the hole levels have been exhausted. One possibil- 
ity is not to carry out a truncation of the hole levels at 
its last iteration, keeping 4p hole states in subsequent it- 
erations. The other is to carry out a truncation at the 
last hole iteration and then keep only p hole states subse- 
quently. The former is more accurate; the latter reduces 
(perhaps significantly) the storage needs for the calcu- 
lation. Considering that storage issues are likely to be 
especially important for subsequent applications of the 
methodology to more complex systems, it is useful to as- 
sess the accuracy lost if we keep only p hole states after 
the hole levels have been exhausted. This is one of the 
issues we will address in the test applications to follow. 

C. Neutrons and Protons 

In the presence of both neutrons and protons, the same 
formalism applies, but with a few minor practical modi- 
fications. 

• In a system of neutrons and protons, there are four 
distinct blocks - neutron particle, proton particle, 
neutron hole and proton hole. 

• In a given iteration, we have a choice of how many 
and which levels to include. One possibility is to 
add a particle level and a hole level of the same 
type, switching between proton levels and neutron 
levels from one iteration to the next. The other 
possibility is to add four levels in an iteration - one 
for each of the four blocks. Both options can readily 
be implemented with the formalism of the previous 
subsection and both permit us to eventually include 
all active levels. We expect that the first approach, 
in which we add a particle and a hole level in each 
iteration, should suffice for heavy nuclei where the 
active neutrons and protons occupy different ma- 
jor shells. For lighter nuclei where they occupy the 
same major shells, it is most likely preferable to 
add all four levels at the same time. 

• The superblock basis in which the hamiltonian is 
diagonalized consists of products of states of all four 



types - neutron particle, neutron hole, proton par- 
ticle and proton hole. Only those states with the 
correct numbers of neutrons and protons and with 
total angular projection zero need be considered. 

• The reduced density matrix for a given block re- 
quires a contraction over the states of the other 
three blocks. 

• We need to consider hamiltonians that include 
single-particle energy terms for neutrons and pro- 
tons as well as the neutron-neutron, proton-proton 
and neutron-proton interactions. The matrix ele- 
ments for all of these terms can be readily obtained 
from the same basic like-particle operators as be- 
fore (at, a^a, a<cft, a) a) a and a^a^aa). 

III. THE MODEL 

The first application of the p-h DMRG, as described 
in the previous section, was to a problem of a very large 
number of particles interacting via a pairing force and 
constrained to a set of equally-spaced doubly-degenerate 
single-particle levels ||. This problem, which is of direct 
relevance to the physics of ultrasmall superconducting 
grains, was solved at half filling for up to 400 levels. De- 
spite the enormity of the model space, especially for a 
very large number of levels, the DMRG method was able 
to very accurately reproduce the ground state of the sys- 
tem, which could be obtained exactly using a method 
pioneered by Richardson jlO) . These results suggest that 
the p-h DMRG method can very accurately describe the 
properties of fermion systems that are dominated by a 
single collective degree of freedom. 

Nuclei, on the other hand, are characterized by sev- 
eral competing collective degrees of freedom, most impor- 
tantly those associated with the pairing and quadrupole 
fields. Thus, before applying the DMRG method with 
confidence to nuclear systems, it is necessary to first 
demonstrate that it continues to work well in the pres- 
ence of such competing collective effects. 

With that in mind, we report here a series of test cal- 
culations for systems that admit pairing plus quadrupole 
correlations on the same footing. The simplest such sys- 
tem is one involving identical nucleons restricted to a 
single level with large degeneracy and interacting via a 
pairing plus quadruple force. Unfortunately, the (0 + ) 
ground state for an even number of particles occupying a 
single level invariably has equal population of all m sub- 
states and thus such a problem does not admit a Fermi 
surface. Since a Fermi surface is an important charac- 
teristic of real Fermi systems, including nuclei, it is im- 
portant to permit one in the test calculations. Thus, in 
the calculations to follow we supplement the pairing plus 
quadrupole interaction of our test model with a single- 
particle energy term that splits the multiply-degenerate 
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large-j oprbit into a series of doubly-degenerate levels. 
The hamiltonian we use takes the form 



H = - X Q ■ Q - flP f P - e E H a] m a 3 



(10) 



The last term splits the levels of the single-j shell into 
a set of equally spaced levels of "oblate" character, with 
the largest \m\ value lowest. Because of the last term, 
the hamiltonian is not in general rotationally invariant 
and thus its eigenstates do not have conserved angular 
momentum. Rotational invariance is of course recovered 
for e = 0. 

In the calculations we will describe, we use a scaled 
version of the quadrupole-quadrupole strength x, defined 
according to 



X 2j + 1 A 



(11) 



where x 1S the usual Q ■ Q coupling strength. With this 
definition, x and g have the same dimensions, which we 
will subsequently ignore. 



In our earlier paper, we presented results for two 
cases, one with x = 1, g = 0, e = 0.1 and one with 
X = 1, g = 0.5, e = 0.1. In both cases, we achieved 
accuracy for the ground state energy of roughly 1 part in 
10 6 with reasonable values of p. 

Obviously the case with g = involves pure 
quadrupole correlations. On the basis of HFB calcu- 
lations that we have now performed, we conclude that 
the g = 0.5 case likewise is dominated by quadrupole 
correlations, with pairing correlations minimal. Since it 
is our desire here to assess the DMRG method in the 
typical nuclear scenario in which different correlation ef- 
fects compete, we will show throughout this section re- 
sults for a further-enhanced pairing interaction, namely 
for x = 1, g = 0.1, e = 0.1. In this case, the HFB calcu- 
lations show a well developed superconducting solution. 
This can be readily seen by noting that for this problem 
the ratio of the average pairing gap to the mean spacing 
between HF levels is 0.39/0.18 = 2.1. 

Figure 2 shows results for the DMRG correlation en- 
ergies, defined as the gain in energy relative to Hartree 
Fock approximation. The exact and HFB values are 



IV. RESULTS FOR 10 PARTICLES IN A J=25/2 
ORBIT 



EcorT 1 = -0.70633 



7^ 



HFB 



-0.20641 



(12) 



A. General remarks 

In this section, we report the results of calculations for 
a system of 10 identical particles in a j = 25/2 orbit. In 
this case, there are 13 active doubly-degenerate orbits, 5 
for holes and 8 for particles. The Fermi surface lies be- 
tween the m = 17/2 hole level and the m = 15/2 particle. 
For this problem, the size of the space that would be re- 
quired for exact diagonalization (assuming that we only 
consider m = states) is 109,583, well within the limits 
of standard diagonalization routines. Thus, for this set of 
calculations, we can readily compare the results of calcu- 
lations based on the DMRG methodology with the cor- 
responding exact results, and thereby obtain important 
insight into the usefulness of the method. We will first 
report calculations for the standard scenario, in which (i) 
we only target the ground state in the DMRG truncation 
strategy, and (ii) we include the full complement of Ap 
hole states in the hole block after all hole levels have been 
treated. We will then present results for scenarios that 
systematically modify these two assumptions. 

B. Results for standard scenario 

Here we present results in which we only target the 
ground state and we include all Ap hole states in the hole 
block after the hole levels have been exhausted. 
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FIG. 2. DMRG correlation energies as a function ofp for a 
system of 10 particles in a j — 25/2 orbit subject to a hamil- 
tonian with x = 1; 9 = 0.1, and e = 0.1. The solid line 
represents an exponential fit to the DMRG results. 
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By a p value of 200, we obtain more than 99% of the 
exact correlation energy, a dramatic improvement over 
HFB approximation which only gets 28%. 

For p = 200 the largest superblock matrix that had to 
be diagonalized had a dimension of 2, 886. For p = 40, 
the maximum dimension was 232; nevertheless, for this 
value of p we still achieve over 95% of the full correlation 
energy. 

Table 1 shows results for the excitation energies of the 
lowest states of this same system. Despite the fact that 
we only targeted the ground state in the density matrix 
phase of these calculations, the agreement is as good for 
the excitation energies of these low- lying states as it is for 
the ground state energy. By p — 200, all three excitation 
energies are reproduced to 1% or better. 

In Figure 3, we show the exact, DMRG and HFB re- 
sults for the occupation numbers of the particle and hole 
levels associated with the ground state solution. The 
DMRG results, which are shown for p = 200, arc in ex- 
cellent agreement with the exact results for all levels with 
appreciable occupation. In contrast, the HFB calculation 
does not obtain sufficient depletion of the Fermi sea. 




From these results, we conclude that the p-h DMRG 
method is able to describe with extreme accuracy the 
low-lying properties of complex many-body systems with 
competing collective features. 



C. Dependence on the number of states targeted 

What is the importance of including excited states in 
the targeting procedure? We address this in Figure 4, 
where we compare the results for the ground state en- 
ergy and the energies of low-lying excited states when 
different sets of states are targeted. The solid triangles 
refer to results when only the ground state is targeted 
(L = 1). The open circles are the results obtained when 
the lowest four states are all targeted in the optimization 
procedure (L — 4). We also include in the figure the 
exact results. 
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FIG. 3. The occupation numbers for the same system as in 
figure 2. The horizontal label 2m is twice the z-projection of 
the angular momentum. The solid circles ( connected to guide 
the eye by the solid line) represent the exact results; the open 
diamonds represent the DMRG results for p = 200; the open 
squares represent the HFB results. 



FIG. 4. The energy of the ground state and of the three 
lowest excited states for the same system as in figure 2. The 
solid triangles are the results when only the ground state is 
targeted in the optimization procedure; the open circles refer 
to calculations when the lowest four states of the system are 
targeted simultaneously. The exact results are represented by 
dashed lines. 
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The results are not unexpected. The ground state is 
better described when it is the only state targeted in 
the optimization procedure. Furthermore, in such a sce- 
nario a reasonable description of low-lying excited states 
is achieved. When several states are targeted, the agree- 
ment gets worse, albeit only marginally, for the ground 
state, but becomes better for the excited states. What 
is most important, however, is that once p becomes large 
enough, there is no discernable difference between the 
two sets of results, both of which are also in excellent 
agreement with the exact eigen-energies. 

D. Dependence on how we treat the hole block after 
all hole levels have been exhausted 

Next we address the effect of maintaining all Ap hole 
states following treatment of the last hole level. This 
was addressed by carrying out the DMRG calculations 
in two ways, namely keeping only p hole states after the 
last hole iteration or keeping Ap. For a value of p = 200, 
we find for the ground state energy a value of —16.7331 
when only p states are retained and a value of —16.7333 
when all Ap are kept. The gain in energy through the im- 
proved treatment of the hole states is just 1 part in 10 6 . 
On the basis of these results, the subsequent calculations 
we describe will be based on the simplifying, but highly 
justified, assumption that we maintain just p hole states 
in those iterations that only add particle levels. 

V. RESULTS FOR 40 PARTICLES IN A J=99/2 
ORBIT 

The excellent quality of the results obtained in the test 
calculations of the previous section has encouraged us to 
treat a much larger system, one for which exact diago- 
nalization is not possible. Following an optimization of 
the DMRG calculational methodology, we are now able 
to treat systems significantly larger than were feasible in 
ref . H . Here we report the largest calculation we have so 
far carried out - for a system of 40 particles occupying a 
j = 99/2 orbit. In this case, the exact calculation would 
involve a hamiltonian matrix of dimension 3.84 x 10 25 , 
obviously much too large to treat without dramatic trun- 
cation. The parameters of the hamiltonian of eq. jl(i| ) 
that were used in this calculation are: \ = 1> 9 — 0-1 an d 
e = 0.2. The calculations targeted only the ground state 
and retained only p hole states after treating the 20 hole 
levels. [The latter is especially important when treating 
extremely large problems such as this one.] The results 
for the ground state correlation energy are shown in fig- 
ure 7. By p = 100, we obtain a ground state correlation 
energy of - 2.8994. An exponential fit to the DMRG re- 
sults, likewise indicated in the figure, give an asymptotic 
result for the ground state energy of — 2.8990 ± 0.0003. 



This is close enough to our DMRG result with p — 100 
to suggest that we have achieved accuracy for the ground 
state correlation energy to 1 part in 10 4 . 

We have also carried out an HFB calculation for this 
system, which confirms that the ground state is super- 
conducting. The correlation energy achieved in this cal- 
culation is — 1.7902, less than 62% of that estimated from 
the exponential fit to the DMRG results. 



-2.80 




I i I i I i I i I 
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FIG. 5. The ground-state energy for a system of 40 identi- 
cal nucleons in a j — 99/2 orbit with the hamiltonian param- 
eters x — 1) 9 — 0-1 and e = 0.2. An exponential fit to these 
results is also plotted. 



VI. CONCLUSION 

In this paper, we describe the recently-developed 
particle-hole Density Matrix Renormalization Group 
method and report its test application to a problem in- 
volving many identical nucleons constrained to a single-j 
shell and subject to a hamiltonian with a pairing plus 
quadrupole interaction and a single-particle term that 
splits the shell into doublets. Earlier work along these 
lines was reported in ref. ||. Here we develop the for- 
malism in greater detail and report calculations that ad- 
dress some important aspects of the method that needed 
to be clarified before the method could be meaningfully 
applied to more realistic nuclear systems. 

In our view, the current calculations confirm and ex- 
pand on the conclusions reported in ref. H regarding 
the use of the method in large-scale shell-model calcula- 
tions. As long as there is a well defined Fermi surface 
in the problem, the method is able to produce extremely 
accurate results for the ground state energy of the sys- 
tem and for the energies of low-lying excited states as 
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well, even in the presence of competing collective cor- 
relation effects. It is also able to accurately reproduce 
other (number-conserving) properties of the system, in- 
cluding for example occupation numbers. Based on an 
optimized computational methodology, we are now able 
to treat much larger problems than in ref. 0], and this 
bodes well for subsequent more-realistic applications of 
the methodology. In all cases we have studied, even those 
very large-scale problems we are now able to handle only 
because of the optimized methodology, we seem to be 
able to obtain accurate results while diagonalizing ma- 
trices of moderate dimensions. Critical to the success of 
the method is the rapid exponential convergence that it 
typically exhibits as a function of the number of particle 
and hole states maintained in a given calculation. 

Clearly the next step is for us to further develop the 
method for use in realistic nuclear systems. This will 
mean including both neutrons and protons, several non- 
generate single-particle levels, and general nuclear hamil- 
tonians. The framework for these extensions is spelled 
out in Section II and is currently in the process of being 
implemented. 
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Table 1: Excitation energies for 10 particles in a 
j = 25/2 level. The Hamiltonian parameters are: x = P 
g = 0.1 and e = 0.1. 



p 


Ei 


E 2 


E 3 


40 


0.89040 


1.26036 


1.49354 


80 


0.88266 


1.18608 


1.41182 


120 


0.88791 


1.18077 


1.41795 


160 


0.88794 


1.17722 


1.41456 


200 


0.88784 


1.17302 


1.41024 


Exact 


0.88578 


1.16245 


1.39915 
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